Anova tables
map(models_q2, ~anova.lme(.x, type = "marginal"))
## $nodule_mass_ratio
## numDF denDF F-value p-value
## (Intercept) 1 35 0.2300 0.6345
## treatment 3 35 0.8484 0.4769
## init_height 1 35 1.2452 0.2721
##
## $nodule_weight
## numDF denDF F-value p-value
## (Intercept) 1 41 3.298 0.0767
## treatment 3 41 4.623 0.0071
## init_height 1 41 1.376 0.2476
##
## $nodule_count
## numDF denDF F-value p-value
## (Intercept) 1 41 27.527 <.0001
## treatment 3 41 11.330 <.0001
## init_height 1 41 5.359 0.0257
##
## $wue_nlme
## numDF denDF F-value p-value
## (Intercept) 1 102 121.59 <.0001
## nfixer 1 6 0.34 0.5798
## treatment 3 102 66.94 <.0001
## init_height 1 102 0.63 0.4307
## nfixer:treatment 3 102 49.70 <.0001
##
## $wue_log
## numDF denDF F-value p-value
## (Intercept) 1 102 36.36 <.0001
## nfixer 1 6 4.09 0.0897
## treatment 3 102 26.88 <.0001
## init_height 1 102 0.40 0.5264
## nfixer:treatment 3 102 4.67 0.0042
##
## $pnue_nlme
## numDF denDF F-value p-value
## (Intercept) 1 102 89.71 <.0001
## nfixer 1 6 0.41 0.5455
## treatment 3 102 16.00 <.0001
## init_height 1 102 7.83 0.0061
## nfixer:treatment 3 102 6.70 0.0004
##
## $pnue_log
## numDF denDF F-value p-value
## (Intercept) 1 102 947.5 <.0001
## nfixer 1 6 0.5 0.4992
## treatment 3 102 3.6 0.0156
## init_height 1 102 2.2 0.1425
## nfixer:treatment 3 102 1.8 0.1541
##
## $n_area_nlme
## numDF denDF F-value p-value
## (Intercept) 1 102 15.472 0.0002
## nfixer 1 6 11.956 0.0135
## treatment 3 102 20.985 <.0001
## init_height 1 102 1.109 0.2947
## nfixer:treatment 3 102 2.026 0.1149
##
## $n_area_log
## numDF denDF F-value p-value
## (Intercept) 1 102 0.606 0.4381
## nfixer 1 6 15.964 0.0072
## treatment 3 102 7.372 0.0002
## init_height 1 102 0.575 0.4501
## nfixer:treatment 3 102 2.385 0.0735
##
## $sla
## numDF denDF F-value p-value
## (Intercept) 1 102 118.25 <.0001
## nfixer 1 6 0.43 0.5342
## treatment 3 102 0.25 0.8618
## init_height 1 102 6.32 0.0135
## nfixer:treatment 3 102 1.80 0.1521
##
## $gs
## numDF denDF F-value p-value
## (Intercept) 1 102 8.874 0.0036
## nfixer 1 6 6.401 0.0447
## treatment 3 102 18.662 <.0001
## init_height 1 102 1.015 0.3161
## nfixer:treatment 3 102 1.698 0.1722
##
## $amax
## numDF denDF F-value p-value
## (Intercept) 1 102 12.364 0.0007
## nfixer 1 6 7.969 0.0302
## treatment 3 102 0.257 0.8562
## init_height 1 102 5.555 0.0203
## nfixer:treatment 3 102 7.384 0.0002
Post-Hoc: Tukey’s test
Maximal photosynthesis
as_tibble(emmeans(models_q2$amax,
pairwise ~ treatment*nfixer,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| no_additions nonfixer - plus_nutrients nonfixer |
-0.1287 |
0.5156 |
102 |
-0.2495 |
1.0000 |
| no_additions nonfixer - plus_water nonfixer |
-0.4310 |
0.5208 |
102 |
-0.8277 |
0.9911 |
| no_additions nonfixer - plus_water_nutrients nonfixer |
-0.0716 |
0.4988 |
102 |
-0.1436 |
1.0000 |
| no_additions nonfixer - no_additions fixer |
-9.0888 |
3.2195 |
6 |
-2.8230 |
0.2419 |
| no_additions nonfixer - plus_nutrients fixer |
-7.0199 |
3.2267 |
6 |
-2.1756 |
0.4646 |
| no_additions nonfixer - plus_water fixer |
-11.4230 |
3.2202 |
6 |
-3.5473 |
0.1120 |
| no_additions nonfixer - plus_water_nutrients fixer |
-9.7834 |
3.2227 |
6 |
-3.0357 |
0.1931 |
| plus_nutrients nonfixer - plus_water nonfixer |
-0.3024 |
0.5379 |
102 |
-0.5621 |
0.9992 |
| plus_nutrients nonfixer - plus_water_nutrients nonfixer |
0.0570 |
0.5175 |
102 |
0.1102 |
1.0000 |
| plus_nutrients nonfixer - no_additions fixer |
-8.9601 |
3.2232 |
6 |
-2.7799 |
0.2532 |
| plus_nutrients nonfixer - plus_nutrients fixer |
-6.8913 |
3.2296 |
6 |
-2.1338 |
0.4829 |
| plus_nutrients nonfixer - plus_water fixer |
-11.2943 |
3.2243 |
6 |
-3.5029 |
0.1174 |
| plus_nutrients nonfixer - plus_water_nutrients fixer |
-9.6548 |
3.2280 |
6 |
-2.9910 |
0.2025 |
| plus_water nonfixer - plus_water_nutrients nonfixer |
0.3594 |
0.5209 |
102 |
0.6900 |
0.9971 |
| plus_water nonfixer - no_additions fixer |
-8.6577 |
3.2230 |
6 |
-2.6863 |
0.2792 |
| plus_water nonfixer - plus_nutrients fixer |
-6.5889 |
3.2302 |
6 |
-2.0398 |
0.5257 |
| plus_water nonfixer - plus_water fixer |
-10.9919 |
3.2235 |
6 |
-3.4099 |
0.1296 |
| plus_water nonfixer - plus_water_nutrients fixer |
-9.3524 |
3.2260 |
6 |
-2.8991 |
0.2232 |
| plus_water_nutrients nonfixer - no_additions fixer |
-9.0171 |
3.2193 |
6 |
-2.8009 |
0.2476 |
| plus_water_nutrients nonfixer - plus_nutrients fixer |
-6.9483 |
3.2267 |
6 |
-2.1534 |
0.4743 |
| plus_water_nutrients nonfixer - plus_water fixer |
-11.3513 |
3.2198 |
6 |
-3.5255 |
0.1146 |
| plus_water_nutrients nonfixer - plus_water_nutrients fixer |
-9.7118 |
3.2220 |
6 |
-3.0142 |
0.1975 |
| no_additions fixer - plus_nutrients fixer |
2.0688 |
0.7107 |
102 |
2.9108 |
0.0810 |
| no_additions fixer - plus_water fixer |
-2.3342 |
0.6752 |
102 |
-3.4573 |
0.0175 |
| no_additions fixer - plus_water_nutrients fixer |
-0.6947 |
0.6816 |
102 |
-1.0191 |
0.9705 |
| plus_nutrients fixer - plus_water fixer |
-4.4030 |
0.7131 |
102 |
-6.1741 |
0.0000 |
| plus_nutrients fixer - plus_water_nutrients fixer |
-2.7635 |
0.7236 |
102 |
-3.8192 |
0.0055 |
| plus_water fixer - plus_water_nutrients fixer |
1.6395 |
0.6758 |
102 |
2.4259 |
0.2403 |
# Treatment effects
emmeans_table_tidy(models_q2$amax,
formula = "treatment|nfixer",
grouping_var = "nfixer")
## [1] "Formula for pairwise comparisons: "
## pairwise ~ treatment | nfixer
## <environment: 0x5f61ade21730>
Water Use Efficiency
as_tibble(emmeans(models_q2$wue_log,
pairwise ~ treatment*nfixer,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| no_additions nonfixer - plus_nutrients nonfixer |
0.0227 |
0.1046 |
102 |
0.2167 |
1.0000 |
| no_additions nonfixer - plus_water nonfixer |
0.6581 |
0.1059 |
102 |
6.2117 |
0.0000 |
| no_additions nonfixer - plus_water_nutrients nonfixer |
0.6950 |
0.1015 |
102 |
6.8499 |
0.0000 |
| no_additions nonfixer - no_additions fixer |
-0.4350 |
0.2152 |
6 |
-2.0215 |
0.5342 |
| no_additions nonfixer - plus_nutrients fixer |
-0.6114 |
0.2197 |
6 |
-2.7829 |
0.2524 |
| no_additions nonfixer - plus_water fixer |
-0.3098 |
0.2154 |
6 |
-1.4384 |
0.8145 |
| no_additions nonfixer - plus_water_nutrients fixer |
-0.2779 |
0.2164 |
6 |
-1.2845 |
0.8769 |
| plus_nutrients nonfixer - plus_water nonfixer |
0.6354 |
0.1092 |
102 |
5.8190 |
0.0000 |
| plus_nutrients nonfixer - plus_water_nutrients nonfixer |
0.6723 |
0.1049 |
102 |
6.4101 |
0.0000 |
| plus_nutrients nonfixer - no_additions fixer |
-0.4577 |
0.2171 |
6 |
-2.1083 |
0.4943 |
| plus_nutrients nonfixer - plus_nutrients fixer |
-0.6340 |
0.2213 |
6 |
-2.8656 |
0.2313 |
| plus_nutrients nonfixer - plus_water fixer |
-0.3325 |
0.2175 |
6 |
-1.5291 |
0.7735 |
| plus_nutrients nonfixer - plus_water_nutrients fixer |
-0.3006 |
0.2189 |
6 |
-1.3736 |
0.8421 |
| plus_water nonfixer - plus_water_nutrients nonfixer |
0.0369 |
0.1059 |
102 |
0.3484 |
1.0000 |
| plus_water nonfixer - no_additions fixer |
-1.0931 |
0.2173 |
6 |
-5.0308 |
0.0254 |
| plus_water nonfixer - plus_nutrients fixer |
-1.2694 |
0.2218 |
6 |
-5.7235 |
0.0136 |
| plus_water nonfixer - plus_water fixer |
-0.9679 |
0.2175 |
6 |
-4.4509 |
0.0443 |
| plus_water nonfixer - plus_water_nutrients fixer |
-0.9360 |
0.2183 |
6 |
-4.2867 |
0.0522 |
| plus_water_nutrients nonfixer - no_additions fixer |
-1.1300 |
0.2151 |
6 |
-5.2535 |
0.0207 |
| plus_water_nutrients nonfixer - plus_nutrients fixer |
-1.3063 |
0.2197 |
6 |
-5.9465 |
0.0113 |
| plus_water_nutrients nonfixer - plus_water fixer |
-1.0048 |
0.2153 |
6 |
-4.6678 |
0.0358 |
| plus_water_nutrients nonfixer - plus_water_nutrients fixer |
-0.9729 |
0.2161 |
6 |
-4.5018 |
0.0421 |
| no_additions fixer - plus_nutrients fixer |
-0.1764 |
0.1445 |
102 |
-1.2206 |
0.9240 |
| no_additions fixer - plus_water fixer |
0.1252 |
0.1372 |
102 |
0.9120 |
0.9843 |
| no_additions fixer - plus_water_nutrients fixer |
0.1571 |
0.1380 |
102 |
1.1379 |
0.9468 |
| plus_nutrients fixer - plus_water fixer |
0.3015 |
0.1448 |
102 |
2.0828 |
0.4330 |
| plus_nutrients fixer - plus_water_nutrients fixer |
0.3334 |
0.1461 |
102 |
2.2824 |
0.3135 |
| plus_water fixer - plus_water_nutrients fixer |
0.0319 |
0.1372 |
102 |
0.2325 |
1.0000 |
# Treatment effects
emmeans_table_tidy(models_q2$wue_log,
formula = "treatment|nfixer",
grouping_var = "nfixer")
## [1] "Formula for pairwise comparisons: "
## pairwise ~ treatment | nfixer
## <environment: 0x5f61aba70700>
as_tibble(emmeans(models_q2$wue_log,
pairwise ~ nfixer,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| nonfixer - fixer |
-0.7525 |
0.189 |
6 |
-3.981 |
0.0073 |
# Treatment effects
as.data.frame(emmeans::emmeans(models_q2$wue_log,
specs = pairwise ~nfixer,
type = "response",
adjust = "tukey")$emmeans) %>%
janitor::clean_names() %>%
dplyr::select(response, everything(),
# Remove variables
-c(df, lower_cl, upper_cl, se)) %>%
# Rename response to emmean, this is done when models is log
dplyr::rename_all(funs(stringr::str_replace_all(., "response", "emmean"))) %>%
# Calculate % difference between control and variable, this assume that
# that first name is the control
dplyr::mutate(difference = ((emmean - first(emmean))),
perc_difference =((emmean - first(emmean) )/first(emmean))*100) %>%
dplyr::mutate_if(is.numeric, round, 3)
## emmean nfixer difference perc_difference
## 1 1.888 nonfixer 0.000 0.0
## 2 4.008 fixer 2.119 112.2
Stomatal Conductance
as_tibble(emmeans(models_q2$gs,
pairwise ~ treatment,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| no_additions - plus_nutrients |
0.0282 |
0.0247 |
102 |
1.1424 |
0.6643 |
| no_additions - plus_water |
-0.1245 |
0.0239 |
102 |
-5.1985 |
0.0000 |
| no_additions - plus_water_nutrients |
-0.1041 |
0.0236 |
102 |
-4.4076 |
0.0002 |
| plus_nutrients - plus_water |
-0.1527 |
0.0251 |
102 |
-6.0899 |
0.0000 |
| plus_nutrients - plus_water_nutrients |
-0.1323 |
0.0248 |
102 |
-5.3291 |
0.0000 |
| plus_water - plus_water_nutrients |
0.0204 |
0.0239 |
102 |
0.8507 |
0.8300 |
# Treatment effects
emmeans_table_tidy(models_q2$gs,
formula = "treatment",
)
## [1] "Formula for pairwise comparisons: "
## pairwise ~ treatment
## <environment: 0x5f61a605f958>
SLA
SLA did not responded to any of the treatments nor varied between nfixers and non-fixing species so no post hoc test was performed.
PNUE
as_tibble(emmeans(models_q2$pnue_log,
pairwise ~ treatment,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| no_additions - plus_nutrients |
0.2100 |
0.0745 |
102 |
2.8193 |
0.0291 |
| no_additions - plus_water |
-0.1149 |
0.0722 |
102 |
-1.5899 |
0.3889 |
| no_additions - plus_water_nutrients |
0.0079 |
0.0715 |
102 |
0.1111 |
0.9995 |
| plus_nutrients - plus_water |
-0.3248 |
0.0758 |
102 |
-4.2844 |
0.0002 |
| plus_nutrients - plus_water_nutrients |
-0.2020 |
0.0755 |
102 |
-2.6757 |
0.0425 |
| plus_water - plus_water_nutrients |
0.1228 |
0.0722 |
102 |
1.6998 |
0.3291 |
# Treatment effects
emmeans_table_tidy(models_q2$pnue_log,
formula = "treatment",
)
## [1] "Formula for pairwise comparisons: "
## pairwise ~ treatment
## <environment: 0x5f61ad643908>
Nitrogen concentration per unit of area
as_tibble(emmeans(models_q2$n_area_log,
pairwise ~ treatment,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| no_additions - plus_nutrients |
-0.1575 |
0.0480 |
102 |
-3.2823 |
0.0076 |
| no_additions - plus_water |
-0.0236 |
0.0465 |
102 |
-0.5074 |
0.9572 |
| no_additions - plus_water_nutrients |
-0.0502 |
0.0461 |
102 |
-1.0895 |
0.6968 |
| plus_nutrients - plus_water |
0.1339 |
0.0489 |
102 |
2.7383 |
0.0361 |
| plus_nutrients - plus_water_nutrients |
0.1073 |
0.0489 |
102 |
2.1925 |
0.1322 |
| plus_water - plus_water_nutrients |
-0.0266 |
0.0465 |
102 |
-0.5722 |
0.9402 |
# Treatment effects
emmeans_table_tidy(models_q2$n_area_log,
formula = "treatment",
)
## [1] "Formula for pairwise comparisons: "
## pairwise ~ treatment
## <environment: 0x5f61b0b25188>
as_tibble(emmeans(models_q2$n_area_log,
pairwise ~ nfixer,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| nonfixer - fixer |
-0.7214 |
0.1866 |
6 |
-3.865 |
0.0083 |
# Treatment effects
as.data.frame(emmeans::emmeans(models_q2$n_area_log,
specs = pairwise ~nfixer,
type = "response",
adjust = "tukey")$emmeans) %>%
janitor::clean_names() %>%
dplyr::select(response, everything(),
# Remove variables
-c(df, lower_cl, upper_cl, se)) %>%
# Rename response to emmean, this is done when models is log
dplyr::rename_all(funs(stringr::str_replace_all(., "response", "emmean"))) %>%
# Calculate % difference between control and variable, this assume that
# that first name is the control
dplyr::mutate(difference = ((emmean - first(emmean))),
perc_difference =((emmean - first(emmean) )/first(emmean))*100) %>%
dplyr::mutate_if(is.numeric, round, 3)
## emmean nfixer difference perc_difference
## 1 1.037 nonfixer 0.000 0.0
## 2 2.133 fixer 1.096 105.7
Nodule weight
as_tibble(emmeans(models_q2$nodule_weight,
pairwise ~ treatment,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| no_additions - plus_nutrients |
-0.0631 |
0.0231 |
41 |
-2.7291 |
0.0443 |
| no_additions - plus_water |
-0.0227 |
0.0129 |
41 |
-1.7523 |
0.3107 |
| no_additions - plus_water_nutrients |
-0.0605 |
0.0198 |
41 |
-3.0512 |
0.0200 |
| plus_nutrients - plus_water |
0.0404 |
0.0239 |
41 |
1.6944 |
0.3397 |
| plus_nutrients - plus_water_nutrients |
0.0026 |
0.0265 |
41 |
0.0978 |
0.9997 |
| plus_water - plus_water_nutrients |
-0.0378 |
0.0205 |
41 |
-1.8480 |
0.2663 |
# Treatment effects
emmeans_table_tidy(models_q2$nodule_weight,
formula = "treatment")
## [1] "Formula for pairwise comparisons: "
## pairwise ~ treatment
## <environment: 0x5f61ad9f9038>
Nodule Count
as_tibble(emmeans(models_q2$nodule_count,
pairwise ~ treatment,
adjust = "tukey"
)$contrast) %>%
mutate(across(2:6, round, 6)) %>%
kable()
| no_additions - plus_nutrients |
0.1686 |
0.0642 |
41 |
2.6250 |
0.0565 |
| no_additions - plus_water |
0.1165 |
0.0689 |
41 |
1.6909 |
0.3415 |
| no_additions - plus_water_nutrients |
-0.1858 |
0.0630 |
41 |
-2.9504 |
0.0258 |
| plus_nutrients - plus_water |
-0.0522 |
0.0648 |
41 |
-0.8059 |
0.8513 |
| plus_nutrients - plus_water_nutrients |
-0.3545 |
0.0640 |
41 |
-5.5420 |
0.0000 |
| plus_water - plus_water_nutrients |
-0.3023 |
0.0683 |
41 |
-4.4277 |
0.0004 |
# Treatment effects
emmeans_table_tidy(models_q2$nodule_count,
formula = "treatment")
## [1] "Formula for pairwise comparisons: "
## pairwise ~ treatment
## <environment: 0x5f61b0a5f530>
Nodule mass ratio
Nodule mass ratio did not responded to any of the treatments applied so no post hoc test was performed.
Boxplots
Boxplots traits
# Step done for getting predictions from models for Q2
data_for_predictions <-
data_for_models %>%
rownames_to_column("id") %>%
# Remove unused variables
dplyr::select(id, spcode, treatment, nfixer, init_height)
string <- c("models_q2")
data_pred_traits <-
# Get models prediction
gather_predictions(data_for_predictions ,
# Return predictions
models_q2$amax,
models_q2$sla,
models_q2$gs,
models_q2$wue_log,
models_q2$pnue_log,
models_q2$n_area_log
) %>%
pivot_wider(names_from = model, values_from = pred) %>%
rename_all(funs(
# rename columns
stringr::str_to_lower(.) %>%
stringr::str_replace(., c(string),"pred_") %>%
# Remove dollar sing
gsub("\\$", "", .)
)) %>%
# Back transform log variables
mutate(pred_wue = exp(pred_wue_log),
pred_n_area = exp(pred_n_area_log),
pred_pnue = exp(pred_pnue_log),
) %>%
# Remove log predictions and init height
dplyr::select(-c(init_height, pred_wue_log, pred_n_area_log))
# Generate plot combinations
vars_q2_interaction <-
crossing(
# Get all numeric variables to plot (all y)
as_tibble(t(combn(dplyr::select(data_pred_traits, where(is.numeric)) %>% names, 1))),
# Select factor variables to plot
x_axis_var = dplyr::select(data_pred_traits, nfixer) %>% names,
group_var = dplyr::select(data_pred_traits, treatment) %>% names) %>%
filter(V1 %in% c('pred_amax', 'pred_wue'))
vars_q2_interaction %>%
# Gererate plots
pmap( ~ boxplot_plot_pmap(data = data_pred_traits,
y = !!sym(..1), x = !!sym(..2),
fill = !!sym(..3)))
## [1] 24.08
## [1] 7.034
## [[1]]

##
## [[2]]

vars_q2_treatment <-
crossing(
# Get all numeric variables to plot (all y)
as_tibble(t(combn(dplyr::select(data_pred_traits, where(is.numeric)) %>% names, 1))),
# Select factor variables to plot
x_axis_var = dplyr::select(data_pred_traits, treatment) %>% names,
group_var = dplyr::select(data_pred_traits, treatment) %>% names) %>%
filter(V1 %in% c('pred_gs', 'pred_n_area', 'pred_pnue' ))
vars_q2_treatment %>%
# Gererate plots
pmap( ~ boxplot_plot_pmap(data = data_pred_traits,
y = !!sym(..1), x = !!sym(..2),
fill = !!sym(..3)))
## [1] 0.359
## [1] 3.399
## [1] 117
## [[1]]

##
## [[2]]

##
## [[3]]

Boxplots Nodules
data_for_nodule_predictions <-
data_nodules_cleaned %>%
rownames_to_column("id") %>%
# Remove unused variables
dplyr::select(id, spcode, treatment, init_height)
string <- c("models_q2")
data_pred_nodules <-
# Get models prediction
gather_predictions(data_for_nodule_predictions,
# Return predictions
models_q2$nodule_count,
models_q2$nodule_mass_ratio,
models_q2$nodule_weight) %>%
pivot_wider(names_from = model, values_from = pred) %>%
rename_all(funs(
# rename columns
stringr::str_to_lower(.) %>%
stringr::str_replace(., c(string),"pred_") %>%
# Remove dollar sing
gsub("\\$", "", .)
)) %>%
# Back transform log variables
mutate(pred_nodule_count = exp(pred_nodule_count),
pred_nodule_mass_ratio = exp(pred_nodule_mass_ratio)) %>%
# Remove log predictions and init height
dplyr::select(-c(init_height))
# Generate plot combinations
vars_q2_nodules <-
crossing(
# Get all numeric variables to plot (all y)
as_tibble(t(combn(dplyr::select(data_pred_nodules, where(is.numeric)) %>% names, 1))),
# Select factor variables to plot
x_axis_var = dplyr::select(data_pred_nodules, treatment) %>% names,
group_var = dplyr::select(data_pred_nodules, treatment) %>% names)
vars_q2_nodules %>%
# Gererate plots
pmap( ~ boxplot_plot_pmap(data = data_pred_nodules,
y = !!sym(..1), x = !!sym(..2),
fill = !!sym(..3)))
## [1] 167.2
## [1] 2.618
## [1] 0.199
## [[1]]

##
## [[2]]

##
## [[3]]
